Evolutionary tradeoff and equilibrium in an aquatic 
predator-prey system 

Laura E. Jones* 
Stephen P. Ellner^ 

February 9, 2008 



Cornell University 
Department of Ecology and Evolutionary Biology 
E145 Corson Hall 
Ithaca, New York 14853-2701 



*lej4@cornelLedu 
^ spe2@cornelLedu 



1 



Abstract 



Due to the conventional distinction between ecological (rapid) and evolutionary (slow) timescales, 
ecological and population models to date have typically ignored the effects of evolution. Yet the 
potential for rapid evolutionary change has been recently established and may be critical to under- 
standing how populations adapt to changing environments. In this paper we examine the relation- 
ship between ecological and evolutionary dynamics, focusing on a well-studied experimental aquatic 
predator-prey system (Fussmann et al. 2000; Shertzer et al. 2002; Yoshida et al. 2003). Major 
properties of predator-prey cycles in this system are determined by ongoing evolutionary dynamics 
in the prey population. Under some conditions, however, the populations tend to apparently sta- 
ble steady-state densities. These are the subject of the present paper. We examine a previously 
developed model for the system, to determine how evolution shapes properties of the equilibria, 
in particular the number and identity of coexisting prey genotypes. We then apply these results 
to explore how evolutionary dynamics can shape the responses of the system to "management": 
externally imposed alterations in conditions. Specifically, we compare the behavior of the system 
including evolutionary dynamics, with predictions that would be made if the potential for rapid 
evolutionary change is neglected. Finally, we posit some simple experiments to verify our predic- 
tion that evolution can have significant qualitative effects on observed population-level responses to 
changing conditions. 



1 Introduction 



A distinction is often made between ecological and evolutionary time scales, the former re- 
ferring to relatively rapid changes in the distribution and abundance of species, the latter to 
more gradual changes in the properties of those species as a result of natural selection. This 
distinction is tacitly embedded in much of the established theory in ecology and evolution. 
Ecological models conventionally assume constant parameters to characterize species and 
their interactions (a recent influential text on community ecology (Morin 1999) docs not 
even list evolution among the "factors influencing interactions among species"). Conversely, 
many evolutionary and population-genetic models assume constant population size on the 
grounds that (relative to the evolutionary time scale) ongoing changes in population size are 
merely high-frequency "noise". When population dynamics are considered, only rare ex- 
treme events (bottlenecks, expansion after colonization, etc.) are assumed to matter. Even 
in theories that explicitly combine ecological and evolutionary dynamics, the separation of 
time-scales is typically assumed. For example, Khibnik and Kondrashov (1997) analyze 
predator-prey coevolution using singular perturbation theory, assuming that ecological dy- 
namics are much faster than evolution. The "canonical equation" of Adaptive Dynamics 
(Marrow et al. 1996, Dieckmann et al. 1997) used in many recent papers (e.g. Dercole 
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et al. 2003, Le Galliard et al. 2003) goes even further, and assumes evolutionary changes 
are so rare that ecological dynamics reach their new asymptotic steady-state or attractor 
before another evolutionary step occurs. 

However, recent years have seen an accumulation of evidence that this presumed separation 
of time scales is often violated (e.g., Thompson 1998, Hairston et al. 1999, Cousyn et 
al. 2001, Reznick and Ghalambor 2001, Grant and Grant 2002, Yoshida et al. 2003): in 
natural and experimental settings, trait evolution and changes in abundance may occur 
on similar time scales. For example, Resniek et al. (1997) observed significant life history 
evolution in guppies over periods of 4-11 years in response to changes in prcdation pressure, 
and estimated that evolution occurred in these populations at rates "up to seven orders of 
magnitude greater than rates inferred from the paleontological record" (Resniek et al. 1997, 
p. 1934). Ashley et al. (2003) review evidence for rapid evolutionary changes in response 
to natural and anthropogenic agents in birds, fishes, mammals, lizards, and plants, with 
changes occuring in some cases within as little as a single year or generation. 

In addition to its theoretical implications, the potential for rapid evolutionary change is 
critical to understanding how populations adapt to changing environments, for example 
in harvested or endangered populations where rapid changes in conditions can result from 
human management or lack thereof (Ashley et al 2003, Stockwell et al. 2003, Zimmer 
2003). A change in extrinsic factors affecting mortality or reproductive success - such as 
a change in harvesting rate or juvenile survival - is necessarily a change in the forces of 
natural selection, and the response is likely to be modified by adaptation. Yet, based on the 
assumed separation of time scales (and despite evidence that commercial fishing has already 
led to significant evolutionary changes in harvested species), plans for species management 
rarely take into account the likelihood that ecological changes will evoke an evolutionary 
response (StockweU et al. 2003, Ashley et al. 2003). 

In this paper, we continue our studies of the interplay between ecological and evolutionary 
dynamics, focused on an experimental aquatic predator-prey system (Pussmann et al. 2000; 
Shertzer et al. 2002; Yoshida et al. 2003). Previous studies (summarized below) have 
shown that major properties of predator-prey oscillations in this system are determined by 

ongoing evolutionary dynamics in the prey population. Under other conditions, however, 
the predator and prey tend to apparently stable steady-state densities. These are the 
subject of the present paper. We study the previously developed model for our system 
under steady-state conditions, to determine how evolution shapes properties of the steady 
states. Specifically, we ask: 

1. Under what conditions do we expect to find single versus multiple genotypes main- 
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tained by selection, and which ones? 

2. To what extent does evolution affect the responses of the system to changes in external 
conditions such as resource availability and extrinsic sources of mortality? 

The idea that evolution can affect ecological dynamics is certainly not new, indeed a sub- 
stantial body of theory (reviewed by Abrams 2000) has developed over several decades to 
predict how the dynamics and stability of predator-prey interactions can be affected by 
trait evolution, and how selective harvesting can affect evolution in fish populations (Law 
and Grey 1989). However, key assumptions of the models have rarely been verified ex- 
perimentally, and the practical difficulties of monitoring long-term changes in population 
abundances and trait distributions in order to test model predictions "has left a rather un- 
fortunate gap between theory and experiment" (Abrams 2000, p. 98). "Although over 40 
years of theory have addressed how evolutionary processes can affect the ecology of predator 
prey interactions, few empirical data have addressed the same issue" (Johnson and Agrawal 
2003). In the absence of an empirical foundation, there has been a proliferation of specu- 
lative models and a corresponding proliferation of alternative predictions: "Evolution can 
stabilize or destabilize interactions . . . When population cycles exist, adaptation may either 
increase or decrease the amplitude of those cycles" (Abrams 2000). In a general model of a 
3-species (resource, consumer, predator) food chain, Abrams and Vos (in press) have shown 
that, as a result of adaptive changes in the consumer, an increase in consumer mortality 
rate could result in either an increase or a decrease in the resource, consumer (prey), and 
predator steady state values. 

Our motivation for adding to this literature is to develop specific predictions for an exper- 
imentally tractable system, where model assumptions and predictions can be tested with 
quantitative rigor, allowing direct feedbacks between theory and experiment as in our pre- 
vious studies with this system (Fussmann et al 2000, Shertzer et al. 2002, Yoshida et 
al. 2003). Experimental tests of theoretical predictions are strongest when predictions 
have been put on record before the experiments are done (Nelson G. Hairston, Sr., per- 
sonal communication). We begin by describing the experimental system and model. We 
then analyze the model using approaches from evolutionary game theory - that is, under 
a given set of experimental conditions we seek to identify either single ESS prey genotypes 
(ESS=evolutionarily stable strategy) that cannot be invaded by any alternative genotype, 
or a non-invasible ESC of coexisting genotypes (ESC=evolutionarily stable combination). 
We then use the results of this analysis to predict how population response to changing 
conditions will be modified by adaptive changes in the prey. 

The model can produce a wide range of dynamical behaviors depending on parameters. A 
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complete analysis of its behavior would be quite lengthy, but also not germane for empirical 
studies. We therefore consider a wide range of values for parameters that are under exper- 
imental control, but for parameters characterizing the organisms we restrict the analysis 
to values near those estimated to hold in our system. Note that experimental conditions 
leading to population cycles require a totally different analysis which will be the subject of 
a later paper. 

1.1 Experimental system 

The experimental system is a predator-prey microcosm with rotifers, Brachionus calyci- 
florus, and their algal prey, Chlorella vulgaris, cultured together in nitrogen-limited, con- 
tinuous flow-through chemostat systems. Brachionus in the wild are facultatively sexual, 
but because sexually produced eggs wash out of the chemostat before offspring hatch, our 
rotifer cultures have evolved to be entirely parthenogenic (Fussmann et al. 2003). The algae 
also reproduce asexually (Pickett-Heaps 1975), so evolutionary change occurs as a result of 
changes in the relative frequency of different algal clones. 

Two parameters of the system can be set experimentally: the concentration of limiting 
nutrient [nitrate] in the inflowing culture medium Nj, and the dilution rate S, the fraction 
of chemostat medium that is replaced each day. A simple ordinary differential equation 
model given below is able to capture the experimentally observed qualitative behavior of 
the system: equilibria at low dilution rates, followed by cycling, followed by equilibria and 
then extinction (Fussmann et al. 2000) as the dilution rate is increased. However, the 
experimental system exhibited longer-period cycles and unique phase relations which did 
not match the short-period cycles and classic quarter-cycle predator-prey phase relations 
predicted by the model. Furthermore, the observed cycles showed extended periods where 
algal densities were high yet rotifer numbers remained low, followed by rapid growth in the 
rotifer populations while algal densities remained roughly unchanged. 

A series of models were devised to explain these observations (Shertzer et al. 2002), each 
including some biologically plausible mechanism: rotifer self-limitation via reduced egg fit- 
ness when food is scarce; changes in algal nutrient composition as a function of nutrient 
availability; changes in algal physiology due to accumulation of toxins released by rotifers, 
and evolved prey defense against rotifer prcdation. Only the model including prey evolution 
was able to account well for the qualitative properties of the observed cycles. In that model, 
algae exposed to rotifer predation pressure evolved a "low palatability" phenotype at some 
cost to their competitive fitness. Subsequent experiments confirmed the existence of an evo- 
lutionary tradeoff between defense against predation and "fitness" (growth rate in absence 
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Table 1. Model Parameters 



Parameter 


Description 


Value 


Reference 


Ni 


Concentration of limiting nutrient in 
supplied medium 


80/i mol N// 


Set 


5 


Chemostat dilution rate 


variable (rf) 


Set 


V 


Chemostat volume 


0.33? 


set 


Xc 


Algal conversion efficiency 


variable 


fitted 


Xb 


Rotifer conversion efficiency 


variable 


fitted 


m 


Rotifer mortality 


0.055/rf 


F2000 


A 


Rotifer senescence rate 


0.4/rf 


F2000 




Minimum algal half-saturation 


4.3 /i mol N // 


F2000 


Kb 


Rotifer half-saturation 


0.835 X 10^ algal cells jl 




Pc 


Maximum algal recruitment rate 




TY 


OJc 


N content in 10^ algal cells 


20.0 11 mol 


F2000 


f-c 


Algal assimilation efficiency 


1 


F2000 


G 


Rotifer maximum consumption rate 


5.0 X lO"'' Ijd 


TY 


a\ 


Shape parameter in algal tradeoff 


variable, ai > 


fitted 


OLl 


Scale parameter in algal tradeoff 


variable, 0:2 > 


fitted 



Set: Adjustable experimental parameter 
F2000: Fussmann et al. 2000 
TY: Yoshida et ai, in prep. 



of rotifers): "grazed" algae grown under constant rotifer predation pressure were smaller, 
competitively inferior, and constituted inferior food (i.e., rotifers grew poorly when fed upon 
these cells) relative to "ungrazed" algae grown in rotifer- free environments (Yoshida et ai, 
submitted). Phenotypic differences between grazed and ungrazed algal lineages were herita- 
ble - persisting in subsequent generations grown under common conditions - demonstrating 
that the algal population evolved in response to grazing pressure. Experiments also verified 
the prediction that predator-prey cycles would exhibit very different qualitative properties 
(shorter period, and different phase relations) in the absence of prey evolution (Yoshida et 
al. 2003). 



2 Description of the Model 

The model is a system of ordinary differential equations describing the predator-prey and 
prey evolutionary dynamics in our microcosms. It is essentially the same as the model used 
by Yoshida et al. (2003), with minor simplifications for the sake of analytic tractability. 
The possibility for genetic variability and evolution in the algal prey population is modeled 
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by explicitly representing the algal population as a finite set of asexually reproducing clones. 
Each clone is characterized by its "food value" to the rotifers, denoted p (for "palatability" ) . 
Food value is defined by its effect on the rotifers: p can be thought of as the conditional 
probability that an algal cell is digested rather than ejected, once it has been ingested. 
The model thus consists of three equations for the limiting nutrient and rotifers, plus q 
equations for a suite of q algal clones. In the following equations, is nitrogen (fi mol 
per liter), Q represents concentration (per liter) of the i^^ algal clone, where i = 1,2, q; 
and R and B are Brachionus, fertile and total population counts in individuals per liter, 
respectively. Fecund rotifers senesce and stop breeding at a rate A; all rotifers are subject 
to fixed mortality m. The parameters Xc Xb are conversions between consumption and 
recruitment rates (additional model parameters are defined in Table I). 



dN 
It 

dCi 
dt 
dR 

H 



d{Nj-N)-^Fc,iiN)Ci (1) 



4 = 1 



XcFc,ii^)Ci - {Ci)B - 5Ci (2) 

XbY,FbM)R-{^ + m + ^)R (3) 

i=l 

^ = XbY.^bSC^)R-{S + m)B (4) 



where 



and 



are functional responses describing algal and rotifer consumption rates, respectively, and 
where pc = 

Equation ^ is derived from the rotifer clearance rate G (the volume of water per unit time 
that an individual filters to obtain food), which in this model is a function of algal food 
value: 

C- G 

Kb + EtidPi' 

That is, algae of lower value result in the rotifers increasing their clearance rate, exactly as 
if the density of food were lower. An alternate model in which rotifer clearance rate does 
not respond to food value was also considered, but could not fit the experimental data on 



FbM 



PcN 
KciPi) + N 

GCm 

Kb + ELi CiPi 



(5) 
(6) 
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A : ai < 1 



B : ai > 1 




0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 

P P 

Figure 1: Plots of clonal food value, p versus half-saturation, Kc{p) for shape parameter m < 1 
(panel A) and ai > 1 (panel B). When the tradeoff curve is concave down (left panel), it is more 
likely that surviving clones, either in equilibrium or in a cycling regime, represent two extremes in 
food value. When the curve is concave up (right panel), intermediate algal types can persist. 

population cycles as well. 



Each algal clone Cj is assigned a food value pi between 1 ("good") and ("bad"), with 
lower Pi resulting in reduced risk of predation. A reduced food value comes, however, at 
the cost of reduced ability to compete for scarce nutrients. This relationship is specified 
by a tradeoff curve [Figure modified from Yoshida et al. 2003. For simplicity, we use a 
two-parameter family to model how Kc varies between Pi = and pi = 1 : 

KM = K, + a2{l-PiT' 

where Kc > is the minimum half-saturation value and ai,a2 > are shape and cost 
parameters, respectively. 

For most parameters we have fairly secure values based on direct experimental measure- 
ments (Table 1). For example, the parameters defining rotifer feeding rate were estimated 
by allowing known numbers of rotifers to graze on algae at known initial densities, with algal 
density measured again after a fixed amount of feeding time (T. Yoshida, unpublished). For 
some parameters, however, we only have indirect estimates obtained by fitting the model 
to population count data - these parameters are described as "fitted" in Table 1. Specifi- 
cally, parameters were chosen to match as closely as possible the amplitude and period of 
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ai < 1 ai > 1 




0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 

deltas deltas 

Figure 2: Numerical bifurcation diagrams for the clonal model. Population densities of algae 
(Chlorella) are shown in grey line; rotifers (Brachionus) in black line. Note the existence of equilibria 
at both low and high dilution rates. Predator-prey cycles occur at intermediate dilution rates. 
Parameter sets used in these calculations were obtained by simultaneous fitting of the clonal model 
[O ~ ®] to amplitude, phase, and phase-lag observations from an intermediate dilution rate regime 
(i.e., S ~ 0.69), and a high dilution rate regime {S = 0.96) during which the system was cycling. 
Left panel. Phase diagram for tradeoff parameter ai < 1: ai = 0.78, 012 = 9.4, and Pmin = 0.10. 
Right panel, Phase diagram for ai > 1: ai = 1.35, a2 = 9.6, and Pmin = 0.10. 



predator and prey cycles, and the phase lag between them, in experimental observations at 
two different dilution rates where the system exhibits cycles {6 = 0.69 and 6 = 0.96), and to 
match the observed lack of cycles at low and high dilution rates. As is often the case, there 
is a range of parameter values more or less equally consistent with the data, and we explore 
how this parameter uncertainty affects our predictions. With the best-fitting parameter 
values, the model (Q-®) produces a bifurcation diagram as a function of dilution rate, 
6, which approximates the behavior our system exhibits in the laboratory (Fussmann et al. 
2000, Yoshida et al. 2003, Supplementary material). "Low-flow" equilibria exist for dilution 
rates 6 ~ 0.05 — 0.50, followed by cycling at intermediate dilution rates, followed by "high 
flow" equilibria at^Ril.O — 1.5 [Figure [2]. 
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3 Analysis 



3.1 Single-clone steady states 

Assume now one dominant algal clone C with food value p. The model (CO)^® then reduces 
to the following four equations, after substituting in the appropriate functional responses 
®-®: 

S(N,-N) 



dt V . / Kc{p) + N 

dC _ PcNC _ GCpB 

'dt~ ^'Kc{p) + N~ (Kb + Cp) ~ 
dR GCpR 

-dt= ^^K^Tc-p-^'^^^'^'' 

dB GCpR ,^ 

-11= + 

After some algebra, the steady state for ((Tj) is found to be: 

N= ^[-7+\/7' + 4iV/^c(p)] 
^ ^ Kb{5 + m + A) 



(7) 



R 
B 



p{XbG - {5 + m + X)) 
5xb{5 + ni)[xc{Ni-N)-C] 
{5 + m + A)2 
8Xb[Xc{Ni -N)-C] 
5 + m + \ 



(8) 



where we define 



7 = Kc(p)-iV/ + ^. (9) 



Note that these only hold for B > If p is too low or 5 too high, the predators will be 
unable to persist, and the steady state nutrient and algal densities are 



iXcpc - S) (10) 
C= XciNi-N) 

3.2 Conditions for ESS and ESC 

We now ask under what conditions a particular clone Or , with an associated food value pr 
might be dominant and non-invasible. Following established usage, a clone C,. is said to 
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be an Evolutionarily Stable Strategy [ESS] if a population consisting of Cr, at steady state, 
cannot be invaded by a rare alternative clone Cj with food value pi. The meaning of "rare" 
here is that the growth rate of a potential "invader" Q is computed on the assumption that 
the population consists entirely of the "resident" type Cr- This quantity is now often called 
the invasion exponent, and is computed as follows: 

WIN dCi XcPcN GpiB 

Ci dt Kc{pi) +N {Kb + CrPr) 

where N, B, Cr are equilibrium conditions set by the resident. Note that X{pi\pi) = 0: 
pitted against itself, a clone neither increases nor decreases. 

We can characterize ESSs by the first and second derivatives of the invasion exponent. To 
be an ESS, an interior trait value p* (i.e., Pmin < p* < Pmax) niust satisfy both the first 
order condition, dX/dpi = and second order condition, d'^X/dpf < at pi = pr = p* ■ A 
trait value satisfying the first order condition will be called a ESS candidate (Ellner and 
Hairston 1994), because the first order condition is necessary but not sufficient. For the 
first order condition we have 

dKPrW) XcPcNaia2{l-pi)^''^-^^ CB 



dpi [K,{pi) + iV]2 (i^, + CrPr 

Setting B* = we see that p* is an ESS candidate if 



(12) 



Taking the derivative of ()12() with respect to the invader's trait value, the second order 
condition is 



-^=i^{2aia2{l-PiT' -{ai-l)[N + K,{pi)]}, (14) 



dpi 
where 

^ XcPciV«ia2(l-p^)^"^~'^ 

provided < < 1. Thus, to be an ESS, an interior candidate p* must also satisfy 

g{p*) = 2aia2{l-p*T' " («i - 1)[^ + K,{p*)] < 0. (15) 

Where no single ESS exists - in particular when a candidate fails to satisfy the second 
order condition - it is often possible to identify a non-invasible ESC (Evolutionarily Stable 
Combination) of coexisting genotypes. An ESC k of coexisting genotypes is characterized 
by the invasion exponent for an introduced rare type pi, 

w I ^ 1 dCi XcPcN CpiB 
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Table 2. Predicted evolutionary outcome as a function of the invasion exponents when oi < 1. 

^{,Pmax\Pmin) ^ ^{Pmax\Pmin) ^ 



^{Pmin\Pmax) ^ 


ESC K {PmimPmax} 


Pmin is ESS 


^{Pmin\Pmax} ^ 


Pmax is ESS 


Both are (local) ESS 



where N ^ B, are the steady state nutrient and rotifer densities set by k and is the 
"effective" algal density at the ESC steady state, 

jeK 

If X{pi\K) < for all Pi ^ K then k is an ESC. A more refined local classification of 
evolutionary stability is possible, also based on derivatives of the invasion exponent (Figure 
1 of Levin and Muller-Landau (2000) summarizes the various stability concepts and their 
relationships), but the ESS and ESC concepts are sufficient for studying our model. 

3.3 Shape parameter ai < 1 

We observe from ()14l) that 

> when ai < 1. (18) 

opf 

Therefore no interior trait p (i.e., Pmin < p < Pmax ) can be an ESS if ai < 1. In addition 
there cannot be an ESC that includes any interior types (see the Appendix), so any ESC 
must consist of the extreme types {pmin,Pmax}- The evolutionary outcome can therefore 
be determined from the two invasion exponents X{pmin\Pmax) and X{pmax\Pmin) [Table 2]. 

Recall that food value p is scaled so that Pmax = !> representing the undefended clone that 
is favored when predators are absent. However, the evolutionary outcome may depend on 
the value of Pmin, the food value corresponding to the highest possible level of defense. We 
consider below a wide range of possible values for Pmim i-e. we regard it as being under 
experimental control rather than a fixed property of the organism. Although the latter 
is literally true, Pmin can be increased temporarily by using a restricted set of founding 
genotypes as in Yoshida et al. (2003). 

Figure IHl shows plots of the invasion exponents X{pmax\p) and X{p\pmax) (panels A and B), 
and the predicted outcome according to the criteria in Table 2 (panels C and D). In the low- 
flow regime (panel A), when p is near both X{pmax\p) and X{p\pmax) are positive. Thus, 
when Pmin is in the range of p values where these inequalities hold, we predict an ESC. 
As p increases, X{p\pmax) remains positive (dashed line) but X{pmax\p) becomes negative 
(solid line). For Pmin in this range of p values, we therefore have X{pmin\Pmax) > and 
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A: Low flow Regime 



B: High flow regime 



CM 

I 



I 





p (clonal food value) 



p (clonal food value) 



C: Low flow regime 



D: High flow regime 



o 



CO 

o 



CM 
O 



O 
O 



ESS is p_min 



:SS is p_max 



1 I I I I r 

0.0 0.2 0.4 0.6 0.8 1.0 
p minimum 





ESS is p_max 




/ ESS is A 


ESC ^ 


/ p_min 1 



"I I I r 

0.0 0.2 0.4 0.6 0.8 
p minimum 



T 
1.0 



Figure 3: Invasion exponents and evolutionary steady states when ai < 1. (A) Plots of the 
invasion exponents A(pmaa;|p) (solid line) and \(j)\pmax) (dashed line) as a function of p for 
low flow rates, b = 0.001 (bold) and 5 = 0.15 (thin). (B) As in panel (A) for high values of 
5,5 = 1.0 (bold) and 5 = 1.7 (thin). Panels (C) and (D) show the predicted evolutionary 
steady states that result from the plots in (A) and (B) using the criteria in Table 2, as 
a function of 5 and pmin- The "corner" in the plot of X{pmax\p) occurs at the minimum 
food value p that allows persistence of the predators. Above this critical value, X{pmax\p) is 
calculated for a system steady state including N, C, and B; below this critical value only 
and C are present in the system at steady state. 
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MPmax\Pmin) < 0, and the ESS is the mmimum existing p (that is, the most defended 
clone). As p increases further the two exponents reverse in sign, so if Pmin lies to the right 
of this change-point Pmax is an ESS. 

The high flow regime is shown in Figure 01 B. For p very near 1 the situation is the same 
as at low flow: X{pmax\p) > and X{p\pmax) < 0, and we again have Pmax the ESS. As 
p decreases X{pmax\p) remains positive for dilution rates at the high end of the high flow 
regime, but becomes immediately negative as p decreases from p near 1 at the low end of the 
regime, i.e., for 5 = 1. X{p\pmax) may either remain negative or become positive, depending 
on the value of 6. When a sign-change occurs, X{p\pmax) is positive and X{pmax\p) negative, 
and we have Pmin as an ESS, or both exponents are positive so the predicted outcome is an 
ESC: a combination of the two extreme types. The location of the sign-change decreases 
with increasing 6, and eventually vanishes. 

The plots shown in Figure 3 are specific to our parameter set, but we show in the Appendix 
that the relevant qualitative properties of these curves are robust to substantial variation 
in parameter values. 

3.4 Shape parameter ai > 1: analysis 

The situation is more complicated if ai > 1, because an internal ESS is then possible. 
We therefore need to analyze (in this Subsection) properties of the first and second order 
conditions for a local ESS. The implications for evolutionary steady states when ai > 1 are 
presented in the next Subsection. 

3.4.1 First order conditions 

We return to the first order expression, ^li{p) (|13|) . and search for ESS candidate types by 
exploring where ()13() holds as a function both of clonal food value p* and dilution rate 5, 
given < p* < 1 and ai > 1. Example calculations of ^i{p) as a function of trait value p, 
dilution rate S, and rotifer conversion efficiency Xb are shown in Figure^ Each zero point 
in a curve corresponds to a ESS candidate associated with a given dilution rate, 5. 

For very low dilution rates, ^li{p) converges to a limiting form shown in bold (Figure^ A, C; 
bold line), for which, in this case, there are three ESS candidates. Depending on the value 
of the cost parameter, 02, the limiting form may either remain below the zero axis (low 02 
values, one ESS candidate), or loop above it (high 02 values, three candidates) as in Figure 
|1]A. As dilution rate increases, the curve dips entirely below the p-axis (Figure |1] A, C; thin 
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A: Low flow regime 



B: High flow regime 




C: Low flow regime D: High flow regime 




Clonal Food value,p Clonal Food value, p 

Figure 4: First order conditions H13|l for the system when ai > 1 (note that cost parameter a2 = 9.8 
can be classified as "high" in value, so conditions for limiting dilution rate (5^0 loop above the zero 
line). A, ili{p) l(T^ . scaled by dilution rate, 6 is shown for the low dilution regime ((5 = 0.1 — 0.5) 
and conversion efficiency Xb = 4000 (thin line). The limiting form for (5 ^ is shown in bold line; 
representative curves for 6 ~ 0.1 and S — 0.5 are shown in thin dashed and solid line, respectively. 

B, fti{p) for high dilution rates S = 1.0 — 1.5, Xb = 4000; curve for 6 = 1.0 is shown in bold; 
representative curves for S — 1.25 and S — 1.5 are shown in thin dashed and solid line, respectively. 

C, As in panel A for S = 0.1 — 0.5, Xb ~ 6000. Again, the limiting form for (5 — > is shown in bold. 

D, As in panel B for (5 = 1.0 — 1.5 and Xb = 6000. Plots for 6 = 1.0 shown in bold. Note that 
in the high dilution rate regime, the ESS candidate values move from left to right as dilution rate 
increases. 
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dashed and solid line), and the two higher-p candidates move together, eventually colliding 
and vanishing. The trait value(s) of ESS candidates are little affected by variations in rotifer 
conversion efficiency Xb in the low dilution regime (Figure |1] A, C). Note that both 6 and Xb 
enter the first order condition, i^i{p), through the steady state values N, B, and C* 

In the high dilution rate regime, there is a single ESS candidate which increases in p-value 
as dilution rate increases (Figure 0] B,D). Lower estimates of rotifer conversion efficiency 
result in higher p-value ESS candidates. For example, at Xb = 4000, an ESS candidate 
exists at p ~ 1 for 8 = 1.0 (Figure I^B, heavy line). However, at Xb = 6000, the candidate 
for the dilution rate 5 = 1.0 sits at p ~ 0.3 (Figure |3]D, heavy line). Higher conversion 
efficiencies result in greater numbers of rotifers, increased selection for defense, and thus a 
lower p-value ESS candidate. 

Most noticeably at higher dilution rates, varying ai also affects the position of the ESS 
candidal. If qi is increased from some initial value ai > 1, candidate p values shift to 
the left. This is sensible when one considers the effect of an increase or decrease in qi on 
the trade-off curve. Increasing ai makes the curve more deeply concave upwards, while 
decreasing ai flattens it. As the curve becomes more deeply concave upwards, one might 
expect surviving clones to migrate in p-value towards the middle of the curve, and as it 
flattens, to move to the extremes - thus, for high p-values, to shift towards the high-p 
endpoint. 

3.4.2 Second order conditions 

Second order conditions (|15|) are strongly affected by variation in several parameters, most 
notably shape parameter ai, rotifer conversion efficiency Xb-, and dilution rate 5. In the 
Appendix we demonstrate that g is an decreasing function of both ai and Xb] a decreasing 
function of 5 in the low dilution range, but an increasing function of 5 in the high dilution 
range. Consequently, increases in either ai or Xb shift the second order curve to the left (at 
a given dilution rate, 5), so that lower p value candidates satisfy the second order condition. 
Conversely, decreasing either oi or Xb will shift the curve to the right, so that only high p 
value candidates, can satisfy the second order condition. 

The behavior of g with increasing 5 implies that ESS candidates must shift to the right 
as a function of 5 in order to satisfy the second order condition (|15() . All else being held 
constant, for ai > 1, only relatively high-p, high-5 candidates thus satisfy both the first 
order ()13() and second order conditions. An internal ESS may exist at high 8 equilibria, 
depending on conversion efficiency Xb^ but not at low 5 equilibria. 
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A: Low flow regime 



B: High flow regime 




n \ \ \ \ \ r n 1 \ \ \ \ r 

0.70 0.75 0.80 0.85 0.90 0.95 1.00 0.70 0.75 0.80 0.85 0.90 0.95 1.00 

clonal food value (p) clonal food value (p) 

Figure 5: Second order conditions ((TH|l for the system when ai > 1. Curves represent second order 
conditions as a function of p* for equilibrium dilution rates S = (0.1 — 0.5; 1.0 — 1.5). A: second 
order conditions for low dilution rate regime 6 = 0.1 — 0.5. Dashed lines show curves for 6 = 0.1, 0.5. 
Heavy bold line shows the limiting form H15|l assumes at very low S (i.e., 6 < O(10~'^)). B: second 
order conditions for the high dilution rate regime (5 = 1.0 — 1.5. In bold line is shown the curve for 
6 = 1.0; dashed line shows curve for (5 = 1.5 

Plotting (fT3|) as a function of p and 6, we see that low-6, low-p candidates invariably fail 
the second order condition, because g{p*) remains positive for these candidates (Figure 121). 
For very low (limiting) 5, there may be high-p candidates which satisfy the second order 
condition; cf. Figure lU panels A, C, bold line. The presence or absence of these candidates 
depends primarily on cost parameter 02, as discussed in the next section. However, as 6 
increases, the ESS candidate shifts to higher p-values, and for these values of p*, g{p*) 
becomes negative. Thus while an internal ESS may exist at high dilution equilibria, there 
is no internal ESS at low dilution rates: a low dilution rate ESS, if it exists, must be an 
endpoint value. 

3.5 Shape parameter ai > 1: Results 

We now use the results above to identify the form of the evolutionary steady states for 
ai > 1. Figure @ summarizes the results. Panels A through D show first and second order 
conditions for this system together on one plot: solid lines indicate ESS candidates, which 
we'll call p*, and the dashed line is a plot of the second order condition g{p*) (|T^ . To the 
left of the dashed line, ESS candidates fail the second order condition; to the right of this 
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line, candidates satisfy this condition. Note that first order conditions for the parameters 
in El A, B are found inHC, D. 

For a low flow regime with high cost parameter value and for -p small, there are ESS 
candidates, p* (solid line, panel A), but they fail the second order condition. Then, for 
Vmin < P* 1 there is thus an ESC of the two extreme types in proportions such that the 
ESC has an average p value approximating the ESS candidate value, p* . For Pmin > P*, 
the first order condition becomes negative and all clones are invasible by Pmin, which forms 
an endpoint ESS. For the high cost- value (q2 = 10) example shown in panel A, there is 
an additional set of high p ESS candidates at very low dilution rate. However, simulations 
indicate that these candidates form only a very local ESS and may be invaded by lower p 
clones. Note that if the cost parameter is lower (e.g., 02 = 4, panel C), these high p, low 6 
ESS candidates do not occur, and only the ESC and low-p endpoint ESS are observed. 

In the high-flow regime, ESS candidates p* generally increase in p-value with increasing di- 
lution rate (panels B and D, solid line). For this parameter set, however, all ESS candidates 
for 5 < 1.35 fail the second order condition. Thus, for Pmin < P*, there is again an ESC 
of 

Pmin and Pmax ) in proportions such that the average p value approximates the candidate 
value p* (panel B). For Pmin > P*, the first order condition becomes negative, so higher 
p clones are invasible by Pmin > P*, and Pmin is again an endpoint ESS. However, where 
ESS candidates do not fail the second order condition, in this case for 6 > 1.35, an internal 
ESS exists (above dotted line, panel B). Recall that as ai is increased, we have shown that 
the second order curve g{p*) shifts to the left. In this example this causes the intersection 
between the first and second order curves (dotted line, panel D), to occur at lower and lower 
dilution rates. Thus more and more of the high dilution rate regime has an internal ESS. 

Recalling that shape parameter ai and assimilation efficiency Xb affect first and second 
order conditions at high dilution rates, we turn to a second example (panels C,D). Here 
both ai and Xb have been increased. Increasing ai shifts ESS candidate values to the left, 
as does increasing Xb- However, this example also features a lower cost parameter (02 ~ 4), 
which steepens the first order line such that all ESS candidate types fail the second order 
condition for these parameters. Thus as above, for Pmin < P* , there is an ESC centered 
about the ESS candidate value, and for Pmin > P*, Pmin forms an endpoint ESS. Since ESS 
candidates fail the second order condition for all dilution rates, an internal ESS never exists. 
The results are otherwise the same as in panel B. 
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A: Low flow regime 



B: High flow regime 
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C: Low flow regime 
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D: High flow regime 
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Figure 6: Evolutionary steady states wlien a\ > 1. ESS candidates are sliown in solid line; 
to the right of the dashed Hne the second order condition |^ < is satisfied, so a candidate 
is an ESS if it hes to the right of the dashed hne. If Pmin lies to the left of the ESS candidate 
line, the evolutionarily stable state is either an ESS (if the candidate satisfies the second 
order condition) or an ESC (if the candidate does not satisfy this condition). Parameters 
for system shown in panels A-B: ai = 1.08, a2 = 9.8, Xb = 5420; panels C-D: ai = 1.19, 
a2 = 3.9, Xb = 6055. 
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4 Simulations 



Since results from the type of ESS analysis we have just completed may only hold locally, 
we tested our analytical conclusions by performing simulations of the system (^3)-®- The 
simulations were written in the R language (version 1.7) and run on a Windows 2000 
platform, using the odesolve package. Tests were performed for a representative sample of 
parameter sets with ai < 1 and ai > 1 and in a range of dilution rates, as follows. 

First, summary diagrams, such as those shown in Figure Ol C,D (for ai < 1) and Figure 
IHl A-D (for ai > 1), were created to show the predicted evolutionary equilibria for each 
parameter set. Then, each diagram was "gridded" coarsely in {pmim and simulations of 
the system in equations Q-® were run to verify the analysis at each grid point. Some 
additional "spot-check" simulations were used to look more closely at predicted transition 
points, e.g. from ESC to ESS. In each case, the system included an initial set of 40 clones, 
distributed evenly in p value along the tradeoff curve determined by the parameter set in 
question. After running a given simulation, we examined time-series trajectories for C, total 
Chlorella cell counts; average palatability p; normalized clonal frequency q as a function of 
time; and mean clonal frequency as a function of p. Each run verified that (1) the system 
was in equilibrium {pmin is a bifurcation parameter, thus equilibrium must be repeatedly 
confirmed as Pmin varies); and (2) the clone or clones identified by analysis were indeed 
dominant. In only one case, noted above, did we find that an ESS identified by our analysis 
was only local in nature, and thus invasible by clones of much lower p value. 



5 Evolution and response to changing conditions 

We can now return to the question posed in the Introduction: how important is an under- 
standing of evolutionary dynamics for making quantitative predictions of how the system 
will respond to "management" (an externally imposed alteration in conditions)? If, as 
the preceding analysis and examples suggest, a change in conditions results in a change 
in community composition, then predictions of how such a community might respond to 
management must not omit the effects of evolution. 

Figure 7 shows predicted responses of the system with and without evolution to changes 
in the dilution rate, which affects both nutrient supply and mortality rates. In panels A,B 
we assume ai < 1, and evolutionary dynamics are as summarized in Figure |2lC. In both 
panels, the dashed line shows the effects on Chlorella density assuming the system responds 
to external change by evolving, while solid line shows the predicted response if effects of 
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Figure 7: Effects of change in dilution rate on steady-state algal density in an evolutionarily 
dynamic system versus a static (non-evolving) system, for ai < 1 (panels A, B) and ai > 1 
(panels C, D). Panels A, B: Pmin = 0.75; solid line shows static response, and dashed line 
shows dynamic response. Parameters and evolutionary dynamics are as in Figure 00. A, 
Initial dilution rate 5 = 0.02, increased to 5 = 0.15. B, Initial dilution rate 5 = 0.15, 
decreased to 6 = 0.02. Panels C, D: parameters and evolutionary dynamics are as in Figure 
IHl D. The dynamic system (dashed line) includes 40 clones; the static system (solid line) 
is a single-clone whose parameters reflect the aggregate properties of the ESC. C, Initial 
dilution rate 6 = 1.1, increased to 6 = 1.4. D, Initial dilution rate 6 = 1.4, decreased to 
6 = 1.1. 
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evolution are ignored. In panel A we set Pmin = 0.75, then slowly increase the dilution rate 
from 6 = 0.02 to 5 = 0.15, thus crossing from "ESS is Pmax^ to "ESS is Pmin^ on Figure 
ISJC. Both systems respond identically initially, but at the "transition" dilution rate {6 = 0.8 
for this choice of parameters), population densities in the evolutionarily dynamic treatment 
suddenly increase to a higher steady-state (dashed line) , while those in the static treatment 
follow the original (solid line) steady-state. In both systems the initial dominant was Pmax, 
but only the evolutionarily dynamic system shifted to the more competitively fit dominant 

type, Pmin- 

In panel B we run the converse experiment: we set Pmin = 0.75, but start at dilution 
rate 5 = 0.15. We then decrease the dilution rate slowly back down to 5 = 0.02, thus 
passing back from "ESS is (Figure |31 C) to "ESS is Pmax" ■ Both systems begin with 

Pmin as dominant and their response is identical, but as dilution rate passes through the 
transition at 5 = 0.08, the evolutionarily dynamic system switches to Pmax as dominant, 
resulting in lower algal densities (dashed line). The evolutionarily static system continues 
to be dominated by Pmin (solid line). So up to a point, the system's response to changing 
conditions can be predicted without regard to the underlying evolutionary dynamics, but 
then a rapid evolutionary response to a small change in conditions causes a discontinuous 
response that would not be predicted when evolutionary dynamics are ignored. 

With ai > 1 (Figure I7p,D) there is an immediate divergence between the actual system 
response, and predictions that ignore the effect of evolution. Again, in both panels, the 
solid line shows the response of an evolutionarily static system, and in dashed line is the 
response of the evolutionarily dynamic, multiple clone system. In panel C we choose an 
initial dilution rate of 5 = 1.10, and slowly increase the dilution rate to 5 = 1.4. At the 
initial dilution rate, the multiple clone system has an ESC comprised of the extremes, Pmin 
and Pmax, in proportions such that their average p value approximates the ESS candidate 
at that dilution rate. Algal densities (dashed line) slowly fall as dilution rate increases, due 
to an increase in average p value of the ESC with increasing dilution rate (higher 6 reduces 
the rotifer density, and thus reduces the selection for unpalatability) . The evolutionarily 
static system has a single clone, whose parameters reflect the aggregate properties of the 
ESC (for example, its value of Kc was set to the nutrient concentration at which the total 
uptake rate of the ESC clone mix reached half its maximum value). Here, higher dilution 
does not evoke the evolutionary response of increased palatability in the prey, so rotifers 
begin to drop out of the system almost immediately, with a consequent steady increase in 
algal densities (solid line). By 6 = 1.35, the rotifers are extinct, and algal densities stabilize. 

We then run our imaginary experiment the other way (Figure [7| D), selecting an initial 
dilution rate of (5 = 1.4, and slowly decreasing the dilution rate back down to 5 = 1.1. 
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The dynamic system responds to a decrease in dilution rate with a slow increase in algal 
densities: this is because the average palatability of clones comprising the ESC is decreasing 

with S, with the expected effect on rotifer abundances. The static system (a single clone 
mimicking the aggregate properties of the ESC at 5 = 1.4) maintains constant palatability as 
conditions change, and again its response is very different from the outcome when evolution 
occurs (note that the system begins to cycle as we pass back through a bifurcation point at 
S = 1.25). 

We now propose a set of experimental checks for our analysis that are robust to parameter 
uncertainties, and in particular to the shape of the tradeoff curve. A general prediction for 
any plausible parameter values is that the evolutionary steady-state at very high dilution 
rates is a single clone (either Pm^a; or an ESS near p„iax), which is replaced by an ESC as 
the dilution rate decreases. As dilution rate decreases further, a Hopf bifurcation occurs 
and the attractor is then a limit cycle. The nature of the cycles (their period and the phase 
lag between predator and prey) depends on which prey clones are present (Yoshida et al. 
2003). 

We therefore predict that if the dilution rate is quickly dropped into the limit-cycle regime, 
the nature of the cycles will (at least initially) depend on the prior dilution rate and evo- 
lutionary steady state. Figure ?? shows this scenario for both ai < 1 (panels A, B) and 
ai > 1 (panels C, D). In panel A (ai < 1), we start with an initial dilution rate of (5 = 1.15. 
The predicted outcome is an ESC of the extreme types, Pmin and Pmax (Figure ?? D). After 
equilibrium is established, we abruptly step down the dilution rate to 6 = 0.70. Given 
the presence of two very different clones, we predict longer cycles with algal and rotifer 
peaks exactly out of phase (Yoshida ct al. 2003). In panel B we start with dilution rate 
5 = 1.65, at which we expect an ESS of Pmax (Figure ?? D). When the dilution rate is 
stepped down to 5 = 0.70, short cycles with classical predator-prey phase relations result, 
as one would expect from a single-clone system (Yoshida et al. 2003). Panels C and D 
show a corresponding step-down experiment for ai > 1. In both cases, the nature of the 
cycles that occur immediately after the step-down reveals the genetic diversity that was 
present before the step-down, and thus lets us test our predictions about the "hidden" prey 
diversity at steady-state. Note, however, that prey evolution continues to occur after the 
step-down (panel D), so only the transient behavior after the step-down is revealing of prior 
evolutionary history. 
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Figure 8: Effects of change in dilution rate on an evolutionarily dynamic system. Panels 
A, B: for ai < 1. Evolutionary dynamics are as in Figure 01 D; parameters are as follows: 
ai = 0.765, a2 = 8.23, Pmin = 0.07, Xb = 6400. Panel A: Initial dilution rate 5 = 1.15 
stepped down to 5 = 0.7 at t = 60; B, Initial dilution rate 5 = 1.65 stepped down to 5 = 0.7 
at t = 60. Panels C, D: for ai > 1. Evolutionary dynamics are as in FigurelHlB; parameters 
ai = 1.087, a2 = 9.79, Pmin = 0.12, Xb = 5419. Rotifers are shown in dashed line and 
Chlorella in bold line. Panel C, Initial dilution rate 5 = 1.2 stepped down to 5 = 0.7 at 
t = 60; D, Initial dilution rate 5 = 1.45 stepped down to 5 = 0.7 at t = 60. 
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6 Discussion 



Our analysis of evolutionary dynamics suggests that small, gradual changes in external 
conditions can precipitate dramatic shifts in community composition or species' densities. 
This can occur for either qualitative type of tradeoff curve (Figure^, through two different 
mechanisms. 

For ai < 1, extreme prey types are favored and an internal ESS cannot exist. Instead, 
the prey population is dominated either by an endpoint ESS (the minimum or maximum 
possible food values), or an ESC consisting of the two extreme types. Under a small change 
in conditions, the evolutionary steady state can jump from one endpoint ESS to another, 
resulting in a discontinuous change in species abundances (Figure [7|^B). 

For ai > 1, an internal ESS may be present at high dilution rates, depending on the shape 
and cost parameters oi and 02, and rotifer conversion efficiency Xb- Where an internal 
ESS does not exist, the population is dominated either by an endpoint ESS composed 
of the minimum food value clone present in the system, Pmim or an ESC of coexisting 
types. Simulations indicate that the ESC consists of the two extreme types PmimPmax in 
proportions such that the average food value of the population approximates that of the 
ESS candidate (which satisfies the first order condition but not the second order condition 
for evolutionary stability). Because the ESC consists of two very different types, the balance 
between them is strongly affected by relatively small changes in conditions. As a result, the 
response of the system (including evolutionary changes in the prey) is very different, not 
only in magnitude but in direction, from what would be predicted if evolutionary changes 
in the prey are ignored. 

Experiments are in progress to determine the shape of the tradeoff curve in our system. 
Circumstantial evidence supports a tradeoff curve where extreme types are favored (ai < 1): 

• The best fits of the model to qualitative properties of the experimental data (in 
particular the transition points between stability and cycles as a function of dilution 
rate, and properties of the cycles) are obtained with ai < 1. 

• The variability in cycle period (roughly 20-40d) observed in experiments can occur in 
the model when ai < 1 (depending on which clones are present), but not for ai > 1 
(Yoshida et al. 2003). 

This raises the possibility that a rapid jump from one extreme prey type to another, as 
predicted by the model, may be experimentally observable in our system. 

Our analysis complements recent work by Abrams and Vos (in press) on a general model 
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for a resource-consumer-predator food chain, with adaptive change in a consumer trait af- 
fecting food consumption and mortality rates. Their goal was to predict how perturbations 
at one trophic level - such as an increase in predator mortality - would propagate through 
the chain and alter other species' abundances, in the presence of consumer adaptation. The 
predictions from purely ecological models that ignore consumer adaptation have sometimes 
been supported, but arc also contradicted by numerous experimental studies (Abrams and 
Vos in press). Abrams and Vos (in press) showed that consumer adaptation broadens the 
range of theoretically possible responses, and so might explain some cases where purely eco- 
logical models failed. For example, without adaptation an incremental increase in consumer 
mortality in their model always entails an increase in resource abundance and a decrease in 
predator abundance; with consumer adaptation the resource abundance may decrease and 
predator density may decrease. 

Abrams and Vos (in press) consider a spectrum of models with different assumptions about 
the composition of the consumer trophic level and the nature of response - e.g., assuming 
a single phenotypically homogenous consumer species responding behaviorally, or a pair of 
consumer species changing in relative abundance - and found that predictions about total 
consumer abundance were robust across a range of models. In our analysis the composition 
of the consumer trophic level is not assumed, but is one of the predicted outcomes from the 
clonal selection dynamics. As a result, we can predict how prey evolution can lead to qiiali- 
tative changes in the genetic composition of the prey population, and how these qualitative 
changes can lead to abrupt population responses to gradual changes in conditions, as well 
as the gradual changes predicted by the Abrams- Vos models. 

Working with well-characterized model ecological systems (e.g., Mueller and Joshi 2000, 
Gushing et al. 2002) makes it feasible to study processes that would be far less tractable in 
the field, and less amenable to rigorous testing of theoretical predictions. Our results lay 
the groundwork for rigorous tests of the longstanding but still contentious hypothesis that 
population management must be "evolutionarily enlightened" (Ashley et al. 2003), rather 
than continuing to take an exclusively ecological perspective (Stockwell et al. 2003). 

A Appendix 

A.l CKi < 1 

To see that we cannot have an ESC including any interior types, consider a coexisting set 
K of two or more genotypes, which has an interior member I and some other member j. By 
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definition we have X{pi\k) = X{pj\K.) = 0, and as a result of the positive second derivative 
(|18|) we conclude that k must be invasible by some other interior genotype near pi. Thus, 
Ki clS defined cannot be an ESC. 

We now analyze the qualitative properties of the invasion exponent functions for ai < 1. 
Consider first X{pmax\p) for low fiow-rate 6 (Figure EK) . For a resident type with p near 
we have B = (the predators cannot persist) so defense against predation has no value 
and any clone with higher p can invade, implying that X{pmax\p) > 0. As the food value 
p of the resident type increases, the predators can then persist (this bifurcation accounts 
for the corner in the plot of X{pmax\p))- We can approximate X{pmax\p) > iov 5 small, by 
inserting the Taylor series expansions of the steady-state values {N, B, C) about 5 = into 
(|TT|) . This gives 

5-^X{p^ax\p) = --ei + 0{p) + 0{6) 
P 

where 9i = NjXciXbG — (m + X)]/[Ki,{m + A)]. The estimated value of 9i from our exper- 
imental data, depending on the value of Xb, is in the range 16.3 — 26.8, causing X(^p.iYiax\p) 
to become negative roughly at p = 9^^ = 0.04 — 0.06. Thus, the qualitative behavior of 
X{pmax\p) at small p is robust to any parameter changes such that ^ 1 continues to hold. 

Similarly, for 6"^ X{p\pmax) near p = we obtain 

5-'X{p\p^ax) = e^^^^^ - {9, -l)p-l + 0{5). 
Kc{p) 

Thus X{p\pmax) > for p small if '^'^^^'j^^ > l/^i, which is true for our estimated parameters 
(the LHS is at least « 0.25, the RHS no more than « 0.06). However, as p increases the 
92 term dominates, and X{p\pmax) becomes negative, as seen in Figure IHK. As (5 ^ both 
scaled invasion exponents approach limiting shapes similar to the curve shown in bold 
(Figure OK). 

For low and high fiow regimes, as p — > 1 the slopes of both invasion exponents become 
infinite (Figure |31^B). With some algebra we can show that this feature depends only on 
the value of ai. In the limit as p ^ Pmax = Ij from ((T^ we have 

dX{pmax\p) , ^ , , dKc 

-I- Constant x — — 

op op 

and 



dp dp 
Since ^ = -aia2(l -ij)"i-^ and ai < 1, then ^^'■P^--\p^ ^ _oo and ^^^^jg"""^ ^ +oo 

as p > Pmax' 
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A.2 ai > 1 



Here we derive the properties of the function g{p) defining the second-order condition, that 
were used in the case ai > 1. As usual these properties depend on parameter values and 
our goal is to identify which parameters or parameter combinations control the relevant 
properties. 

(a) g{p*) is a decreasing function of cti for all cti > 

Substituting the steady state expressions into (fT^ . we rewrite 5 as a function of ai. 
Taking a partial derivative of this expression with respect to ai yields 



dgjp*) _ K,{p)' 
dai 2 



(3ai -1) - (ai -1)- ^ + 2^^ 



[^ + ^c_-a2(l-pr] 



where 7 is defined in equation Q, and K'^{p) = ^^^^^^ = a2ln{l — p){l — p)°''^ . 

Note that K'^{p) is negative for any < p < 1. The expression in the first of the two square 
brackets is always positive and remains bounded between ~ 1.4 — 2 as Nj is either increased 
from its present value by up to ten times or reduced to zero, and Xb assumes the range of 
values we obtained in our optimized parameter sets {xb = 4000 — 6500). The expression in 
the second square bracket has a lower bound of ~ 4.5, and is also always positive. We thus 
have: 

dai L ^ J L ^ i 

and g{p*) is always a decreasing function of ai for any dilution rate 5 within our experi- 
mental range, and for biologically reasonable choices of the experimental parameter Nj and 
the fitted parameter Xb- 

(b) g{p*) is a decreasing function of Xb 

Substituting the steady state expressions (jSJ into (fT3)) . we rewrite 5 as a function of Xb- 
Taking a partial derivative of this expression with respect to Xb yields 

dXb 2 y^-/^ + miKc{p*) 

where 7 is defined in equation @. .^(7) thus has the same sign as 7', which may be written 
as follows: 

, ^ 97 ^ PcKhG{5 + m + A) 
^ dXb ^6p[xbG-{6 + m + X)]^^- 

The expression within the curly brackets is always positive, thus -^(7) < and thus g{p) is 
a decreasing function of Xb for ai > 1. 
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(c) g{p*) is a decreasing function of S for 6 small, and an increasing function of 5 
when 6 is large for ai > 1 

Again substituting the steady state expressions ^ into (|15|) , we can rewrite g' as a function 
of 6. Taking a partial derivative of this expression with respect to 6 yields 

where F{'y) is the expression on the right-hand side of equation (|19() and 7 is defined in 
equation @. As above, F{'y) has the same sign as 7', which may be written as follows: 

, _dj_ _ PcKb ^ + m + - XbGjm + A) 
^ ~ 86 ~ p* ^ [6{xbG - m - A) - 6^^ 

This implies that F^j) < at low 6, and -^(7) > at high 6, switching sign at dilution rate 

^crit = -(m + A) + s/xbGim + A). 
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